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Abstract 


This paper presents a direct approach for the determination of aerodynamic performance 
characteristics of horizontal axis wind turbines. Based on Glauert’s solution of an ideal wind- 
mill along with an exact trigonometric function method, analytical closed form equations are 
derived and given for preliminary determination of the optimum chord and twist distributions. 
The variation of the angle of attack of the relative wind along blade span is then obtained 
directly from a unique equation for a known rotor size and refined blade geometry. A case 
study including the analysis of an existing turbine model is given and results are discussed 
and compared with those obtained by other investigators. It is shown that the approach used 
in this study is efficient and saves much of the computational time as compared with the 
commonly used iterative procedures. © 2001 Elsevier Science Ltd. All rights reserved. 
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Nomenclature 
a Axial induction factor 
a’ Angular induction factor 
Ap 
C Blade chord 
Ci Normal force coefficient 
Ci Tangential force coefficient 
Cp Drag coefficient 
CL Lift coefficient 
Cp Power coefficient 
Co Torque coefficient 
Cr Thrust coefficient 
F Tip-loss factor 
Np Number of blades 
P Extracted power 
Qaro Aerodynamic torque 
r Local blade radius 
R Rotor radius 
Taro Aerodynamic thrust 
Vo Mean wind velocity at hub height 
V; Resultant wind velocity 
Vg Net wind velocity 
a Angle of attack 
Og Airfoil setting angle 
0, Initial built-in twist 
O; Setting angle at blade root (Blade pitch) 
A Tip-speed ratio 
A, Local speed ratio 
p Air density 
(oy Local solidity 
@ Inflow angle 
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1. Introduction 


The success of wind power as an alternate source of energy is a direct function 
of the economics of designing and manufacturing the wind machine. One of the 
most important aspects of the design is the development of performance prediction 
methods. It is the purpose of this study to present analytical and direct procedures 
for evaluating performance characteristics of horizontal-axis wind turbines. Glauert 
[1,2] originated the basic aerodynamic analysis concepts of airscrew propellers and 
windmills. He applied first the momentum and energy relationships for simple axial 
flow and then considered the effects of flow rotation after passing through the rotor 
as well as the secondary flows near the tip and hub regions. Wilson et al. [3,4] 
extended Glauert’s work and presented a step by step procedure for calculating per- 
formance characteristics of wind turbines. Analysis was based on a two-dimensional 
blade-element strip theory and iterative solutions were obtained for the axial and 
rotational induction factors. Later on, Wilson [5] outlined a brief review of the aero- 
dynamics of horizontal-axis wind turbines. The performance limits were presented, 
and a short note was given on the applicability of using the vortex-flow model. 

This paper focuses on the improvement of the methods of performance analysis. 
The suggested approach is divided into two convenient stages, namely; (a) the evalu- 
ation of a preliminary blade geometry, i.e. the chord and twist distributions, and (b) 
the calculations of the various aerodynamic coefficients for a refined blade geometry 
which can be practical from the manufacturing and cost points of view. In the first 
stage a closed-form cubic equation in the axial induction factor is developed accord- 
ing to Glauert’s optimum solution of an ideal windmill [1]. This is then solved 
analytically by an exact trigonometric function method yielding the optimum chord 
and twist distributions. The optimum working operating condition of the airfoil sec- 
tion is considered which leads to the maximum lift-to-drag ratio. Based on strip 
theory equations [4], the second stage derives a unique equation in the angle of 
attack of the resultant wind velocity and the various aerodynamic coefficients are 
then determined systematically. The method is investigated by testing an existing 
turbine model and the results are compared with those obtained by other investi- 
gators. It is shown that the proposed approach simplifies and reduces much of the 
computational efforts as needed by the traditional iterative methods. 


2. Performance-optimized wind turbine 


A wind turbine system design should be based on the cost-effective production 
of energy, which depends very much on the site wind characteristics and machine 
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performance characteristics. Maximization of the annual energy production may be 
attained by maximization of the rotor power coefficient. This depends, for a pre- 
scribed tip-speed ratio, on the type of airfoil, chord and twist distributions and num- 
ber of blades. 

Several theoretical studies have been published for the determination of a perform- 
ance-optimized wind turbine. Glauert [1] initiated the calculation of the optimum 
wind turbine by making the power integral equation stationary. He considered an 
ideal actuator disk model and obtained the optimum variations of the axial and 
rotational induction factors (a and a’), as well as the inflow angle, , as follows 
(refer to eqs. 2.10, 2.11 and 2.14, page 328 in Ref. [1]): 


13a 
" “4al oe 
a'(1+anA2=a(1—a) (1b) 
d-a) 
tan °F 4a’) (Ic) 


where A,=Qr/V,, is the local speed ratio at the rth station along the blade, Q is the 
rotational speed, and V,, is the net wind velocity. Fig. 1 shows a front view of a 
rotating blade and the velocity triangle of an arbitrary airfoil section. In deriving 
Eqs. (la), (1b) and (1c), secondary effects of drag and tip losses were neglected. If 
one further neglects rotation behind the rotor (i.e. a=1/3 and a’=0), the inflow angle 
may be determined from the relation: 


g=tan“'(2/3/,) (2a) 


which may be used only for 2,>1. 
Glauert presented the solutions of Eqs. (1a), (1b) and (1c) in a tabulated form (see 


(a) (b) 


n 
dr 
af 


Fig. |. Airfoil-blade geometry. (a) Front view of rotating blade; (b) Velocity triangle. 
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Table 34, Page 329 of Ref. [1]) by first assigning a value for the factor a, calculate 
the corresponding value of a’, and then the local speed ratio A, It is a major aim 
of the present study to calculate the factors a and a’ directly from analytical 
expressions for a specified value of A,. Wilson et al. [4] performed a local optimiz- 
ation analysis by maximizing the power output at each radial station along the blade. 
The axial induction factor was varied until the power contribution became stationary 
(see page 59 of Ref. [4]). Rohrbach and Worobel [6] investigated the effect of blade 
number and section lift-to-drag ratio on the maximum turbine performance. Their 
results have been found to yield slightly lower maximum performance than that 
found in Ref. [4]. 

Optimum power coefficient was also investigated by Nathan [7] who derived an 
approximate relationship between the inflow angle (@) and the local speed ratio (A,). 
The relation was given by the following Sth order polynomial: 


$=57.51—35.56A,+10.61A2—1.58643+0.11444—0.00313A5 (2b) 


where @ is measured in degrees. Nathan’s equation was obtained for a lift-to-drag 
ratio ranging from 28.6 to 66.6 and the effects of secondary flows in the tip and hub 
regions were not included in the analysis. 

Fig. 2 shows the variation of the optimum inflow angle, @, with the local speed 
ratio as calculated by the various methods. It is remarked that the results obtained 
by Wilson et al. (see Table 3.2.3, page 74 of Ref.[4]) are identical to those determined 
from Glauert’s Eqs. (la, 1b) and (1c). Deviations can only be significant at stations 


Glauert (Œq). 


Equation (2.8). 
Nathan(Eq.2ib). 


Wilson (4). | 


Inflow angle, degs. 


Local speed ratio, A, 


Fig. 2. Variation of optimum inflow angle with speed ratio. 
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very close to the blade tip region. Comparisons between the results obtained by 
Glauert [1] and Nathan [7] indicates that the major discrepancies occur near the 
blade root portion, while no significant difference is observed along the rest of the 
blade. More examination of Fig. 2 indicates that the optimum distribution of the 
inflow angle can be adequately determined from Glauert’s solution. The next section 
will be devoted to the development of a closed-form analytic expression for the 
inflow angle as determined from Glauert’s optimum conditions. 


2.1. Analytical optimum Glauert’s solution (an exact trigonometric function 
method) 


Combining Eqs. (la) and (1b), we get the following cubic polynomial in the axial 
induction factor a: 


16a?—24a? +3(3—-A?)a—(1 -A?) =0 (3) 


Introducing the transformation a=b cosO+1/2, and comparing with the trigono- 
metric identity 4cos*@—3cos@=cos30, we obtain two values for the parameters b and 
6 as follows: 


br=+A/2 (4a) 
g= oos (EA; (4b) 
A,=/1+42 (4c) 


It can be shown that the root which ensure a positive value for the rotational 
induction factor must be given by the relation 


a =U —A,(cos@*—cos@)] (5) 


The factor a’ and the corresponding optimum value of the inflow angle @ can then 
be determined automatically from Eqs. (la) and (1c), respectively. 


2.2. Optimum chord and twist distributions (preliminary blade design) 


Glauert defined the theoretical aerodynamic efficiency of a blade element, 7, as 
the ratio of the power absorbed by the element to the useful work done by the thrust 
on it (refer to page 214, eq. 2.3 of Ref. [1]). This leads to the relation: 


_ I-tan ycot @ 


a l+tan ytan @ (Sa) 
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tan y= C (6b) 
L 


It can be seen that the efficiency increases as the drag-to-lift ratio decreases. Fig. 
3(a) shows the polar diagram of a typical airfoil section. It is shown that tany is 
minimum for the angle of attack, œ, corresponding to the point where the straight 


(a) 


CL 


lift coefficient 


drag coefficient Cp 


Stall due to 
flow separation 


i eo SR 


= 


Ja lift slope 


Fig. 3. Aerodynamic characteristics of a typical airfoil section. (a) Polar diagram; (b) C,—a@ curve. 
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line from the origin becomes tangent to the curve. For this particular value of a, 
the efficiency reaches a maximum value and, hence, the turbine performance. For 
unstalled condition, the relation between Cp and C, may be approximated by a 
second-degree polynomial in the form 


Cp=Cpot KiC tK C? (7a) 


where Cpo is the drag coefficient at zero lift, and K, and K, are constant coefficients 
to be determined from the prescribed data of the airfoil section such as those 
presented in Ref. [8]. 

Therefore the optimum condition can be given by 


CL= y CpolK2 (7b) 


Cp=2Cpo tK CL (7c) 


The corresponding minimum value of tany is equal to K,+2VK-Cpo. The optimum 
angle of attack can then be determined from the relation 
CL 
=+ 
aC daa” (8) 
where œ, is the zero-lift angle of attack and 0C,/dq@ is the C,—a@ curve slope (see 
Fig. 3(b)). 
Referring to Fig. 1, the optimum twist angle can be calculated from the relation 
On=o-a (9) 
and the optimum distribution of blade chord, c, for given number of blades (Np), is 
determined from 


_ Ngc 
OT ir (10) 


a 


where ø, is termed as the local solidity ratio which can be calculated from Eq. (15) 
for optimum condition (see next section). It is to be noticed that, in Eq. (10), the 
dimensionless chord c and radial distance r are defined by the notations c—c/2R and 
r——r/R, respectively. This means, for example, that the dimensionless airfoil chord 
is equal to its dimensional value divided by twice the rotor radius (refer to Table 1). 


Table 1 
Definition of dimensionless quantities 


Quantity Notation Non-dimensionalization 
Radial distance of blade element from hub center r r—r/R 

Radial distance of element from blade root x x—x/L 

Blade root offset from hub center Ty Tyt—lry/R 

Blade length L L—L/R 

Blade chord Ç C—C/2R 


Wind Velocity (Vo Vas Va V—V/OR 
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3. Aerodynamic performance analysis via a direct method 
3.1. Basic equations 


Wind turbine performance is determined in terms of aerodynamic thrust, torque 
and power. The present analysis utilizes the basic equations of strip theory (see eq 
2.15, p. 329, Ref. [1]. Also see eqs 2.4.5—6, p. 23-24 and eq 2.6, p. 39 of Ref. [4]), 
where the axial and rotational induction factors at any radial station, r, along the 
blade are given by 


a Oo, = 1 (11a) 
l-a 4F sin “ 4F sin’ a 
— +1 
0,0, 
i GG, 1 
a sq! (1 1b) 
l+a' 4F sin @ cos @ 4F sing cos 1 


0,C, 


The size and rotational speed of the rotor as well as the distribution of the chord 
and twist are all assumed to be preassigned parameters. The blade geometry may 
be determined from refined optimum distributions as obtained in the previous Sec- 
tion 2.2. 

Ca and C, are known as the normal and tangential aerodynamic coefficients, 
respectively, and can be determined from (see Fig. 1) 


Ca =C, coso + Cp sind (12a) 
C,=C, sing — Cp coso (12b) 


F denotes Prandtle’s tip-loss factor defined by the relations [1,4]: 


pe He” 13 
= (e7) (13a) 
_ Ng (l-r) 

J= 2 rsing (Lab) 


According to Wilson [4], the tip-loss factor is introduced herein as a first order 
contribution, which is referred to as the standard method of analysis. Dividing Eq. 
(11b) by Eq. (11a) and substituting in Eq. (1c), yields to 
Ci a' 

a 


A, (14) 


c7 
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Combining Eqs. (lla) and (11b) with Eq. (14) gives the following transcen- 
dental equation: 


AF sing(cosd—A,sin ¢)—0,(A,C, + C,) =0 (15) 


which represents a unique relation that can be used directly to determine the angle 
of attack and, hence, the various aerodynamic coefficients at the specified blade 
station. Any suitable root-finding interactive-computer routine, such as that based 
upon Newton-Raphson or Bisection method, can be used to find the required sol- 
utions of the above equation. 


3.2. Calculation of aerodynamic thrust, torque and power 


Having determined the variation of the angle of attack and the corresponding 
aerodynamic coefficients (C1, Cp, Cn, C,) along the blade span, the total thrust, torque 
and power developed by the N,-bladed rotor can be computed by summing up the 
contributions of the individual blades and integrating over one complete revolution. 
Firstly, we calculate the dimensionless velocity of the relative wind from the relation 
(refer to Fig. 1 and Table 1). 


Ve= Val ay tar +a')’] (16) 


Secondly, the thrust, torque and power coefficients are determined from 


1 


2NpL | V2CC,dx 
Thrust — Coefficient Cray (17a) 
1 
anat | ViCCr dx 
Torque—Coefficient Co=— va (17b) 
Power— Coefficient Cp=4Co (17c) 


where A=OR/V, is the tip-speed ratio of the rotating blade (refer to Appendix A for 
the details of deriving Eqs. (17a), (17b) and (17c)). Finally, the aerodynamic thrust, 
torque and power can be determined from the relations: 


1 
Thrust Tyero=5PVeACr (18a) 
1 
Torque Qaero = 3P VApRCQ ( 1 8b) 


1 
Extracted Power P=5P VADC = QQ aero (18c) 
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where Ap is the effective disk area swept by the rotating blades and R is the rotor 
radius. 


4. Application (a case study) 


As a case study, the proposed method of analysis will be implemented on an 
experimental, 100 Kw horizontal axis wind turbine; namely, the ERDA-NASA 
MOD-0 model. The machine operates at a design rotational speed of 40 rpm while 
producing its rated power at 8 m/s wind velocity. The rotor has two blades with 
38.1 meters diameter located downwind from the tower at 30 meters height above 
the ground. The blades are preconed at 7° and have NACA 230xx-series airfoils. 
More about the configuration description and technical data can be found in Ref. [9]. 


4.1. Optimum chord and twist distributions 
Based on the analytic procedure outlined in Sections 2.1 and 2.2, the optimum 


chord and twist distributions of the turbine model described above have been determ- 
ined and presented in Table 2. To compare the results with those obtained by Wilson 


Table 2 
Optimum chord and twist distributions for MOD-0 turbine (C,=0.9, ~@=8°, Cp=0) 


Station Optimum Chord (meter) Optimum Twist (degrees) 
r Wilson [4] Present Wilson [4] Present 
0.05 3.4653 3.4673 34.255 34.3289 
0.10 3.5643 3.5749 21.995 22.0489 
0.15 3.0264 3.0384 14.460 14.5050 
0.20 2.5210 2.5321 9.710 9.7490 
0.25 2.1281 2.1378 6.534 6.5679 
0.30 1.8288 1.8374 4.290 4.3192 
0.35 1.5978 1.6054 2.630 2.6560 
0.40 1.4161 1.4226 1.357 1.3804 
0.45 1.2698 1.2757 +0.351 0.3732 
0.50 1.1500 1.1553 -0.462 -0.4413 
0.55 1.0506 1.0549 -1.134 -1.1130 
0.60 0.9665 0.9697 -1.699 -1.6760 
0.65 0.8952 0.8961 -2.182 -2.1546 
0.70 0.8330 0.8309 -2.613 -2.5663 
0.75 0.7794 0.7708 -2.978 -2.9241 
0.80 0.7312 0.7115 -3.325 -3.2380 
0.85 0.6867 0.6460 -3.661 -3.5155 
0.90 0.6422 0.5610 -4.023 -3.7625 
0.95 0.5919 0.4265 -4.494 -3.9839 
1.00 0.0000 0.0000 - = 


Power Coeff Cp 0.547 0.5312 
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et al. [4], the effect of drag is ignored, i.e. Cp is set equal to zero in Eqs. (12a) and 
(12b). All sections composing the blade are assumed to have NACA 23018 profile 
operating near the optimum working condition where C,=0.9 and a@=8°. The constant 
coefficients of Eq. (7a) have been determined from the data presented in Ref. [8] 
for Reynold’s number equal to 3.0x10°. They are: 


Cpo=0.0067985, K,=—0.00105613, and K,=0.0068344. 


Other parameters used in generating the results are: 


C.—a@ Curve Slope: 0C,/da=5.73 per radian. 
Zero-lift angle of attack: @=—1.1° 
Max. lift coefficient: Cimax=1.4 


Moreover, the present analysis neglects the effect of wind shear and tower shadow 
which means that the net wind velocity, V,, at the various blade stations can be 
taken equal to the mean wind velocity at hub height, V.. 

Examining the results presented in Table 2, it is seen that the values obtained by 
Wilson et al. [4] are too close to those calculated by the present method of analysis. 
Deviations can only be significant in the neighborhood of blade tip portion. It is also 
noticed that the difference in calculating the power coefficient (Eqs. (17b) and (17c)) 
is about 2.8%. The actual design case with C,/Cp=79 has also been investigated 
where the power coefficient was found to be 0.4657 which is less than the value of 
0.473 given in Ref. [4] by about 1.5%. Such discrepancies may be attributed to the 
fact that Wilson procedure contained a second-order contribution from the tip-loss 
factor. However, experimental results for wind turbines are not sufficiently accurate 
to ascertain the validity of including the 2nd order term. 


4.2. Performance of a trimmed rotor (constant power condition) 


For a preassigned blade configuration and rotor speed the power output of a wind 
turbine is basically a function of the wind velocity and pitch angle at blade root. 
Wind turbines are usually designed to generate, within their operating range, a con- 
stant amount of power at a specified rotor speed. Therefore, given the mean wind 
velocity at hub height, V., the problem is to determine the proper value of the pitch 
angle, 0, , at root required to maintain a specified constant power output (trimmed 
rotor). Ref. 10 presented a numerical procedure for the determination of the desired 
angle and the computational results as applied to the MOD-0 machine were plotted 
and given in Fig. 4. The total setting angle, 0g, of an arbitrary blade section at station 
x is given by (refer to Fig. 1): 


OB, Vo) =0 Vo) —0,(4) (19a) 


where 0, is the initial built-in twist of the blade. The following exponential function 
is interpolated and found to be a best-fit representation of most of the data given in 
Ref. [9]: 
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Fig. 4. Pitch angle and aerodynamic coefficients of a trimmed rotor. 
0(x)=35.0138[1.0—e°26%] (in degrees) (19b) 


The refined blade chord distribution of the actual MOD-0 design is defined by the 
following dimensionless expressions: 


C(x)=0.036 0.0$x<0.2 


(19c) 
C(x)=0.042-0.03x 0.2Sx=1.0 


Using the developed transcendental Eq. (15), Fig. 5 depicts the variation of the angle 
of attack, œ, along the blade span for wind speeds of 8, 10, 12 and 16 m/s. All fall 
within the operating range of the MOD-O turbine model. In each case the pitch angle, 
O., is adjusted by the pitch-control mechanism so that the rotor develops its rated 
shaft power of about 133 kw (i.e. trimmed rotor). It is observed that œ is well- 
behaved and continuous function in the radial distance x along the blade span. Con- 
vergence to a unique solution has been achieved in a minimum time of computation. 
The present method of analysis eliminates the possibility of having multiple sol- 
utions, which may occur when using the traditional iterative procedure [4]. It is to 
be noticed that a feasible solution of Eq. (15) must ensure positive values of the 
induction factors. The axial induction factor should not exceed, in any case, a value 
of 0.5. Fig. 4 also shows the variation of the calculated aerodynamic coefficients 
with wind velocity. 
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Fig. 5. Variation of angle of attack along blade span for different wind velocities. Trimmed rotor con- 
dition (constant power). 


5. Conclusions 


A new formulation of performance-prediction equations of horizontal-axis wind 
turbines has been presented. These equations were successfully applied to an existing 
machine; namely, the ERDA NASA MOD-0, producing 100 kw electrical power. 
The optimum aerodynamic blade geometry as well as the trimmed-rotor solutions 
were obtained and investigated in detail. Comparisons with the traditional methods 
of analysis have shown that the present approach is computationally efficient and 
ensure continuity of the angle of attack along blade span. This, in turn, guarantees 
convergence and uniqueness of the attained solutions from the newly developed tran- 
scendental equation. Furthermore, it has been indicated that the optimum blade con- 
figuration can be adequately determined from an exact trigonometric function 
method, which is based on Glauert’s optimum conditions. In conclusion, the pro- 
posed method of analysis eliminates much of the numerical efforts as required by 
other iterative procedures. Finally, problems that remain are the inclusion of tower 
shadow and wind shear effects as well as the influence of the choice of the type of 
airfoil section on the overall aerodynamic efficiency. These problems will be hope- 
fully reported in future studies. 
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Appendix A 


It is convenient to deal with dimensionless quantities so that the analysis can be 
valid for arbitrary blade configurations and rotor sizes. The various non-dimensional 
quantities are defined in Table 1. The same symbols that define the actual parameters 
will be reused to define their corresponding non-dimensional quantities in order to 
avoid having more subscripts and symbols in the manuscript. For example, the 
notation V,—V/QOR, as defined in Table 1, means that the dimensionless resultant 
velocity is equal to its dimensional value divided by the tangential velocity at the 
blade tip. 

Applying blade-element theory, the thrust coefficient can be calculated from: 


R 


1 
| sPVICCNpdr 


T, ny 
an = i aero i ( Al ) 
SPVA PVR -r 


Substituting for r=ry+x (dr=dx) and nondimensionalizing the various parameters 
of the integrand 


1 


(ia (Soran 


C=" (A2) 
K a y (2) 
OR) 7 R 


Using the notations defined in Table 1, one must get the same form of Eq. 
(17a), i.e. 


1 
INL | ECC 


Cr= (A3) 


0 
mVo(1—rii) 


Similar derivations shall lead to Eq. (17b). Eq. (17c) is derivable from the defi- 
nition of the power coefficient. 
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